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Related Applications 

This application is a continuation of US Provisional Application Serial Number 
60/213,189 filed June 21, 2000. That disclosure is incorporated in this document by reference. 

Technical Field of the Invention 

The present invention is generally related to resource-planning methods used by 
manufacturing companies and other organizations. More specifically, the present invention 
includes a strategic resource planning method that accounts for uncertainties inherent in the 
forecasting process. 

Background of the Invention 

hi the field of product manufacturing, the term assemble-to-order (ATO) refers to systems 
where products are not produced until the demand for the products becomes known. 
Manufacture-to-stock (MTS) systems, on the other hand, base production on component 
availability. Unlike ATO manufacturing, MTS systems often produce products before the 
demand for those products has been determined. 

Manufacturing companies typically prefer to use ATO manufacturing operations 
whenever possible. ATO allows manufacturers to minimize their component inventories. This 



reduces the amount of capital invested in components and reduces the risk that components will 
loose value or become obsolete before they are transformed and sold as products. 

In ATO manufacturing operations, it is paramount for decision-makers to be able to 
accurately forecast production on a period-by-period basis. The absence of accurate forecasting 
means that makes the estimation of important quantities such as expected revenues and costs 
difficult and in some cases impossible. 

Unfortunately, accurately forecasting production can be problematic. Product demand is 
necessarily an important factor in estimations of this type. At the same time, product demand is 
not, by itself, a sufficient basis for forecasting. This follows because production can (and often is) 
constrained by the availability of the components needed to assemble the products. As a result, 
there is always a chance that the demand for a product will not be met within a particular period. 
Typically, this means that customers are either turned away or that orders are delayed - the 
product is "back ordered." 

To accurately forecast production, manufactures must be able to compute expected 
production from stochastic demand data, component consumption data and component levels. In 
many cases, this computation will involve a plurality of product and component kinds, 
potentially numbering in the thousands or more. Figure 1 describes this problem for a simplified 
case where there are two products PI and P2 and two components CI and C2. 

In Figure 1, the production constraints imposed by the components CI and C2 are shown 
as two dashed lines labeled CI and C2. These lines divide the space of possible production 
values for the pair (PI, P2) into feasible and infeasible regions. The feasible region is the portion 
of the positive quadrant under the lines CI and C2. 

The probability distribution of demand is illustrated by the concentric ellipses. Each 
ellipse is an iso-probability curve (i.e. the line connects (PI, P2) points of equal probability). It is 
clear that there is a nonzero chance that the demand for the product pair may fall in the infeasible 
region. 

If demand falls in the infeasible region, actual production will be less than demand. 
Hence mean production will be less than mean demand. There is a need for systems that can 
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accurately calculate the differences between mean demand and mean production. This is 
particularly true for manufacturing operations that involve large numbers of products or large 
numbers of components. It is also particularly true where the markets for products or the market 
for components are volatile. 

5 Summary of the Invention 

An embodiment of the present invention includes a method for forecasting the mean 
production (the expected production) for a target planning period. To begin this method, a user, 
or planner chooses one or more products for which the expected production is desired. 

The user then enters data describing each selected product. The data entered for each 
£}j product includes data describing the demand for that product as well as data describing the 
^ components required for each product. Inter-product dependencies are also entered. 

:\ The expected production for each selected product for the target planning-period is 

-=-: expressed as a sum of multidimensional integrals involving the data entered in the previous two 

- steps. Once formulated, the integrals are evaluated. The result of this computation is then 

1-5 presented to the user. 

Stated differently, the present invention includes a method for the estimation of mean 
production for assemble-to-order manufacturing operations, the method comprising the steps of: 
specifying a product to be analyzed; entering data describing the components required to produce 
the specified product; formulating a sum of multidimensional integrals corresponding to the 
2 0 estimation of mean production for the specified product; and evaluating the sum of 
multidimensional integrals. 

Other aspects and advantages of the present invention will become apparent from the 
following descriptions and accompanying drawings. 

Brief Description of the Drawings 

2 5 For a more complete understanding of the present invention and for further features and 

advantages, reference is now made to the following description taken in conjunction with the 
accompanying drawings, in which: 
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Figure 1 is a graph showing production as it relates to demand in both feasible and 
infeasible regions. 

Figure 2 is a block diagram of an Internet-like network shown as a representative 
environment for deployment of the present invention. 

5 Figure 3 is a block diagram of a computer system as used within the network of Figure 1 . 

Figure 4 is a flowchart showing the steps associated with an embodiment of the mean 
production forecasting method of the present invention. 

Figure 5 is a flowchart showing the sub-steps associated with the step of entering data 
about products and components as used within the method of Figure 4. 

1 6 Figure 6 is a flowchart showing the sub-steps associated with the step of generating a sum 

.'.I of integrals as used within the method of Figure 4. 

-= Figure 7 is a flowchart showing the sub-steps associated with the step of selecting a 

\ production policy as used within the method of Figure 4. 

Detailed Description of the Preferred Embodiments 

1 5 The preferred embodiments of the present invention and their advantages are best 

understood by referring to Figures 1 through 4 of the drawings. Like numerals are used for like 
and corresponding parts of the various drawings. 

Definitions 

Component plan: a list of quantities for each component, representing a company's 
2 0 component order for a given planning period. 

Scenario: a set of assumptions about products and components. A scenario includes 

product parameters, component parameters, component consumption, component 
interactions, and an allocation policy. 
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Environment 



In Figure 2, a computer network 200 is shown as a representative environment for an 
embodiment of the present invention. Computer network 200 is intended to be representative of 
the complete spectrum of computer network types including Internet and Internet-like networks. 
Computer network 200 includes a number of computer systems, of which computer system 202a 
through 202f are representative. Computer systems 202 are intended to be representative of the 
wide range of large and small computer and computer-like devices that are used in computer 
networks of all types. Computer systems 202 are specifically intended to include non-traditional 
computing devices such as personal digital assistants and web-enabled cellular telephones. 

Figure 3 shows a representative implementation for computer systems 202. Structurally, 
each computer system 202 includes a processor, or processors 300, and a memory 302. Processor 
300 can be selected from a wide range of commercially available or custom types. An input 
device 304 and an output device 306 are connected to processor 300 and memory 302. Input 
device 304 and output device 306 represent all types of I/O devices such as disk drives, 
keyboards, modems, network adapters, printers and displays. Each computer system 202 may 
also includes a disk drive 310 of any suitable disk drive type (equivalently, disk drive 310 maybe 
any non- volatile mass storage system such as "flash" memory). 

Overview of Method and Apparatus for Resource Plan Analysis Under Uncertainty 

As shown in Figure 4, an embodiment of the present invention includes a Method 400 for 
forecasting the mean production (the expected production) for a target planning period. Method 
400 begins with step 402 where a user, or planner, chooses one or more products for which the 
expected production is desired. 

In step 404 the user enters data describing each product selected in step 402. The data 
entered for each product includes data describing the demand for that product as well as data 
describing the components required for each product. Inter-product dependencies are also 
entered. 

In step 406 the expected production for the target planning-period is expressed as a sum 
of multidimensional integrals involving the data entered in the previous two steps. 



In step 408 the sum of integrals formulated in the previous step is evaluated. This 
evaluation may be performed using a range of different methods including Monte Carlo 
simulation and quadratures. It is preferable, however to use the method disclosed in the 
copending commonly owned application " Method And Apparatus For Multivariate Allocation Of 
Resources ." This method is particularly preferable in cases when Method 400 is being performed 
for large numbers of products or components. In these cases, other evaluation methods may be 
too slow or unstable. 

In step 410 the result of the computation of step 408 is presented to the user. 
Product Selection 

In step 402, the user identifies the m products that will be analyzed. 
Product and Component Data Specification 

Step 404, where the user specifies the required data for the mean production calculation, 
can be further subdivided into the sequence of steps shown in Figure 5. The initial portion of this 
sequence is a loop that includes steps 404-2 through 404-10. In this loop, data is input for each of 
the m products that have been selected for processing by Method 400. Data is also input for each 
of the n components that are required (either directly or indirectly) to produce the products that 
have been selected for processing by Method 400 

In the loop of steps 404-2 through 404-10 the following information is entered: 

1) The type of components required to produce each product, 

1) The number of each component required to produce each product (A), 

1) The allocation of each component (d), 

1) The mean demand for each product 

1) The standard deviation of the demand for each product, 

1 ) The correlation between each product and any other product. 
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Within this loop, each data item can be input manually by the user (see step 402-6) or 
retrieved from a user database (see step 402-8). 

In the loop of steps 404-2 through 404-10, an m x n matrix A is created. The matrix A 
includes one column a ; for each of the m products selected for analysis. Each column 
a, contains the bill of materials required to produce one product. Each of the n elements of a 
given a, specifies the quantity of a particular component required to produce the corresponding 
product. 

The loop of steps 404-2 through 404-10 also generates a vector d. The vector d includes a 
total of n elements, one for each component. Each element within d specifies the available 
quantity of a corresponding component. 

The loop of steps 404-2 through 404-10 also generates a vector d. The vector d includes a 
total of n elements, one for each component. Each element within d specifies the available 
quantity of a corresponding component. 

The loop of steps 404-2 through 404-10 also generates a vector \i . The vector \i includes 
a total of m elements, one for each product. Each element within p specifies the mean demand 
for a corresponding product. 

Once the data required for each product and component has been input, processing 
continues at step 404-12. In step 404-12, the product correlations (item 6 above) and demand 
standard deviations (item 4 above) input in the preceding loop are combined to create a 
covariance matrix 2. 

Formulation of Expected Production as a Sum of Multidimensional Integrals 

In Step 406, the covariance matrix is used to formulate a sum of multidimensional 
integrals. This step (i.e., creation of the multidimensional integrals) can be further subdivided 
into the sequence of steps shown in Figure 6. 

In step 406-2, a feasible region (denoted Q) is defined. The feasible region Q is defined to 
include all points where the demand for a product can be met with the current levels of the 
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components needed to manufacture the product. To determine Q, the matrix A and vector d 
generated in step 404 are used. Q, is defined as the space of those q vectors, for which the 
inequality d - Aq > 0 holds, where q is a vector variable denoting the production of the 
products. 

In step 406-4 a production policy is specified for the feasible region Q. The production 
policy for a region is a function that relates production to demand. Within Q, component 
availability d equals or exceeds demand. As a result, production is equivalent to demand within 
the feasible region €1. For this reason, the production policy for Q is set to q(x) = x reflecting 
the fact that production is capable of mirroring demand within Q. 

hi step 406-6, an infeasible region (denoted Q ) is determined. The infeasible region Q is 
the space of all vectors, where the demand for a product exceeds the supply of components 
needed to manufacture the product. Q includes all vectors, or equivalently points, which are not 
included in Q. Q is defined as the complement of Q as determined in step 406-2. In general, Q 
and Q are separated by a hyperplane in the n dimensional space 91* . 

In step 406-8 a production policy for the infeasible region Q is formulated. Unlike 
production in the feasible region Q, production within the infeasible region Q does not mirror 
demand. This is because component demand exceeds component supply within Q . In step 406-8 
a production policy is specified to determine how production will proceed in light of component 
supply inadequacies. The production policy is a function q(x) that maps demand to production in 
fi. 

A range of different production policies maybe used in step 406-8. For the particular 
implementation being described, a production policy known as the uniform production policy is 
used. The uniform production policy relates production q to demand x within the infeasible 
region Q using the function: q(x) = x - a(x) u . In this function, u is a production policy 
vector that the user may optionally specify. The production policy vector u shows the amounts by 
which demand for products has to be decreased until the demand is feasible. «(x) is defined as: 
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, v max [: 



a ; 



a, 



x 



- di 



u 



This formulation of a(x) ensures that, for any x e Q and a given u, the production policy 
a(x) is the highest level such that q(x) < x and d - Aq > 0 . The values of the 
elements of the constraint matrix A can be calculated from the components of the policy vector u 



In step 406-10 the expected or mean production qis expressed as the following sum of 
multidimensional integrals: 



: where f(x, p, Z) is the multivariate normal density function with mean \i and covariance 

1 G = _ matrix E . 

':= Evaluation of Multidimensional Integrals 

: = In step 408 the sum of the integrals, which were formulated in the previous step, is 

evaluated to compute the expected or mean production q . To compute q , a(x) is rewritten 
using the Cholesky decomposition of 2. As previously noted, £ is the covariance matrix 
1 5 originally computed using the product demand standard distributions and inter-product 

correlations (see step 404-12). The Cholesky decomposition of the covariance matrix £ can be 
written as TT'=£, where T is a lower triangular matrix. A benefit of the Cholesky decomposition 
is that the transformed random variables z have zero mean and 7 n , the identity matrix as their 
covariance matrix. With this transformation of variables x may be rewritten as 

2 0 x = Tz + ji , and the feasible region D also gets transformed into Q'. 

The Case of Single Component Constraints 

In this section we demonstrate the essentials of the technique by first describing the case, 
when there is a constraint only for a single component. In this case the feasible region Q is 



5 and a(x). 
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defined by the constraint equation d -a \i-b z = 0, where b = T'a, and d is a scalar number, 
describing the single constraint. 

After the Cholesky transformation the next step is to perform an orthogonal rotation R, 
which transforms the vector b to point along a unit vector of the new system, say ei: Rb = |b| ei, 
where |b| is the magnitude of vector b. After the rotation the new variables will be denoted by w: 
w = Rz. Because of the orthonormality of R the covariance matrix of w will still be the identity 
matrix. An advantage of this rotation is that the hyperplane, defining the feasible region Q is now 
perpendicular to only one axis, defined by e u and it is parallel to all the other axis of the rotated 
coordinate system. Therefore the n dimensional integral, confined by a generic hyperplane, is 
transformed into an integral, which is confined along a single axis, and is unconfined in the other 
n-1 directions. The unconfined integrals are straightforward to execute and we are left with a 
singlevariate integral, which is much easier to perform. In mathematical form the integral of 
some function h(x) over the feasible region can be written as: 



Here h"(w) is the generic function h(x) after the Cholesky transformation T and the 
rotation R, w„.i = (w 2 , . . ., w n ), and is the (n-1) dimensional space of real numbers. Finally 
U = (d- aji) / | b | , identifying the intersection of the hyperplane and the ei axis. 

In the special case of h(x)= d - ux , the twice transformed function takes the form : 

/z"(w) = <i-u p-RT'u w . 

With this the integral I can be rewritten as 

I = (d-u-\i)3[-oo,U]-s l 3, [-<», I/] • 

Here 



u 




u 



3[-oo, U] = J /(w^ ;0, I n _ x )dw „_! J f(w 1 ;0,1) dw t 
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and 

u 

3 x [-°°,C/]= J f(^„-.i;O r I H _ l )dw n _ l J wJiw^O,!) dw x ' 

-co 

finally 

To summarize: as a result of the Cholesky transformation and the rotation R, the original 
I multivariate integral over a confined subspace has been reduced to a univariate definite integral, 
which can be executed often analytically, and always numerically. 

The Case of Multicomponent Constraints 

We now consider the case of having several constraints for the availability of 
components. These constraints can be written as a set of constraint equations 

d - Aq > 0 

where d is ap dimensional vector, and A is &p times n dimensional matrix. Here p is less 
than or equal to n. The i-th row of the matrix A corresponds to the constraint vector aj similarly 
to the case of the single component constraint. The first step of the method, the Cholesky 
transformation, can be performed in the present case as well, transforming to the new variables z. 
The constraint vectors ai are transformed into the vectors bi by this transformation. 

These equations once again define the feasible region, which is now bounded by a set of 
hyperplanes. However in a generic case the hyperplanes are not orthogonal to each other. For this 
reason in general it is hard to separate the multivariate integrals into univariate integrals. 

An approximate evaluation of the integral / can be constructed the following way. Define 
a new r-dimensional orthogonal basis u defined by a set of r hyperplanes, such that each of the 
original m hyperplanes is mapped to a new hyperplane. In general r is less than or equal to m and 
in practice r is often 1. The index mapping from the original m-dimensional basis to the new r- 
dimensional set is represented by a such that cr{i) is the index of new hyperplane associated 
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with the original hyperplane i.. The new basis may be defined in various ways with the objective 
being to closely approximate the originalfeasible region. In practice we take the new basis to be 
a single orthanormally transformed hyperplane. 

In general the dimension of the space of these new u orthogonal constraint vectors can be 
5 smaller than those of the original b constraint vectors. Therefore some of the transformed 
hyperplanes, which are now orthogonal to only one axis, will be parallel to each other, hi this 
case the integration along the say j-th axis should go up to Cj, the minimum of the integrational 
limits along these axis. In terms of the original variables: 



where r < m. 

Now we will employ the above derived relations for the case, when we have to integrate a 
linear function of the variable A(x)=di-ajX . The twice-transformed function takes the form: 



min 




With this finally the integral I is given as: 
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A"(w) = d, - a f ji - RT'a, • w 



Introducing 8 = 1^3, 



With this the integral I can be rewritten as 



7 = (d i -a rl i)n 3Ho, Ci ]-£ ^3J-co fCj ] f]3[-oo, Ci ]. 
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With all these preparations, we are now ready to evaluate the integral of the production 

policy. 

We start with performing the same Cholesky decomposition, which allows a(x) to be 
reduced as: 

,/ \ /™ \ max fb,- • z + a ; • u - d ; ~] 

a'(z) = a(7i + |t) = . r , -2 '- ^ '-\ 

i e \m\ |_ a ; • u J 

where again b ( . = T'a { + fi . 

Next the same orthornormal rotation is performed as before, using the matrix R. Then 
the m vectors b ; are again approximated by r orthonormal vectors u y such that 

b ; = jb, | u CT (,) . This approximation identifies mutually orthogonal hyperplanes, which will 
allow for a more straightforward evaluation of the integral /. As before, a is a mapping from the 
index set [m] into [r]. Substituting for the b ; and defining the rotation matrix R such that 
Ru a ^ = e a ^ where e, are unit vectors yields: 




max 




a, - u 



and y. = 



a. 



- d i 



u 



allows a"(yv) to be rewritten as: 




We are now ready to compute the value of the expected production. An important part of 
the calculation is to evaluate the following integral 



/* = J«(x) f{x, v, S) dx . 

n 

Using the above preparatory formulae I* reduces to 

1=1 -00 

where we used that the fact that the f multivariate normal density function is normalized 
to one, thus in the directions where the integration is unconstrained, the integrals can be readily 
performed to yield unity. Hence only those integrals remain, where a hyperplane introduces a Cj 
upper bound for the integral: 



mm 
j e<r _1 (0 



Since «"(w) is a linear function of its argument, once again the preparatory formulae can 
be applied to yield 



> max / 



Using the notation 



M, 



rj) 



the following describes the integrational limit U ik : 
mm [M,. - Yj 
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In the final analysis the original n dimensional multivariate integral with an arbitrary 
covariance matrix, constrained by a set of p arbitrary hyperplanes, has been reduced to r 
independent univariate definite integrals. Note that when r=l, a~\i) is the entire set [m]. 

As can be shown, the original integral is an NP hard problem, requiring exponentially 
long time for evaluation, whereas the just-derived univariate integrals can be evaluated in a 
straightforward manner by any known numerical technique, including even the Simpson 
procedure. 

The total expected production is then given by 

q = Jx f(x, 11, E) dx + jq(x) f{x, ji, E) dx 

q n 

q = J x f(x, \i, E) dx - u Ja(x) f(x, \i, E) dx 

n+n n 

q = Jx f(x, ji, E) dx - u (a(x) f(x, n, E) dx 

n+n n+n 

+ u Ja(x) f(x, p., E) dx 
n" 

Here one recognizes that the integration over the Q + Q. space is an integration over the 
total space, not constrained by the hyperplanes, and as such, can be readily calculated. Finally the 
last term is exactly what has been calculated above, thus the above result for I* can be directly 
employed here, yielding the desired expected production. 

The Local u-Policy 

The uniform policy prescribes a demand-production mapping that is manifestly irrational 
in some cases. In particular, the problem of "intercomponent" effects arises where a component 
that gates production for one product will diminish production over all. This arises because the 
previous policies compute a single a for all products. 
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On an event by event bases, an "iterated policy" represents rational demand-production 
mappings, but analytic formulas are not available for the expected production over these policies. 
A compromise is to compute a separate a for each product. This we refer to as the local u-Policy. 

As before, for product i, let D (z) denote the components in the BOM of i. Also, as 
5 before, for any component j we can define the subspaces Q = jx : dj - a y • x > oj . We now 
define the subspaces Q 1 for each product i as 

:i If x e Q' then component allocation is sufficient to meet product demand for the i th 

: product, even though it may not be sufficient to meet demand for other products. 

lQn We make the important observation that in computing expected values for functions 

h (x) that can be separated into a sum, 
= h(x) = h l (x 1 ) + ... + h n (x n ) 

.-« the appropriate subspace of integration for the feasible production for each function hj is Q}. 

1 5 For any subspace Q, let Eq denote the expectation over support Q. Thus, the expectation 

ofh can be expressed as 

E* = |KM*) + E ff M*(x))| 



where q(x) is any production policy that maps each x e Q' into q t (x) e . 

20 

We define the local ^-policy q as 



x t if x e Q ( - 

-a ; .(*, otherwise 



2 5 where 
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a,(x,d) 



We now show that with the preceding definition, for any k, we have that d k - a k - q > 0 . 
Substituting q we get 



5 But observe that a kl = 0 if k <£ D(z') , and therefore 
max a. • x - d. 



proving that d R - a k ■ q > 0 . 

1Q-; Computation of the results in this case is identical to the method described in the previous 

section except that a separate integration is done for each element of the expected production 
vector. For each element i, the index set of constraints [m] is replaced by the reduced set D (i) 
which includes only those components in the BOM for product i. (i.e., the index of rows in the i- 
= th column of the connect matrix A for which the entry is nonzero). 

1 5 Although particular embodiments of the present invention have been shown and 

described, it will be obvious to those skilled in the art that changes and modifications may be 
made without departing from the present invention in its broader aspects, and therefore, the 
appended claims are to encompass within their scope all such changes and modifications that fall 
within the true scope of the present invention. 
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1 Expectation of General Functions Under a Multivariate Normal 
Distribution 



We introduce some proprietary results on the properties of of the multivariate normal distribution 
over arbitrary half-spaces. Let /(x; fi, S) denote the multivariate normal density function of x e Jr.™ 
with mean vector fx and covariance matrix £. Let Ct denote the half space fi = {x:o£ — a ■ x > 0} 
and let /i(x) denote an arbitrary multivariate function of x. We show how to obtain an expression 
for the following general integrals: 



where Cl = {x : d - a • x < 0}. 

The solution to ly and 7 2 can be obtained by performing two transformations. We first perform 
a Cholesky decomposition that transforms the multivariate normal distribution into a multivariate 
normal with zero mean and identity matrix for the covariance matrix. The Cholesky decomposi- 
tion amounts to a linear transformation of the space U n , and therefore transforms Cl into another 
halfspace Cl'. After the Cholesky decompositions, we rotate the space so that the hyperplane defin- 
ing the halfspace CI' lies along one of the coordinate axis in the Cholesky decomposition space. 
We exploit the fact that the zero-mean, identity matrix covariance matrix multivariate normal dis- 
tribution is invariant under SO(n) group transformations, where SO(n) denotes the Lie group of 
orthogonal rotations in 5R U . 

1.1 Cholesky Decomposition 

For a positive definite matrix £, the Cholesky decomposition defines a unique lower triangular 
matrix T such that TV = £ [Tong, p 184] . We note that x has distribution AT([i, E) . Let T denote the 
matrix of the Cholesky decomposition of E. It follows from the linear transformation properties of 
multivariate normal distributions that the random variable z = T~ l (x— ji) has distribution M(0, I n ) 
[Tong, p 32]. With that transformation, the defining hyperplane for the subspace fi becomes 




(1.1) 



id the dual integral 




(1.2) 



d - a - (Tz + (m) = 0, 



(1.3) 



where b = T'a. The transformed hyperplane defines the new region $7' = {z : d - a • p, - b • z > 0}. 
The transformation also defines a new function h'(z) = h(Tz + p). Thus, we get 

l x = y ft'(z)/(z;0 ) 7 n ) £ fa (1.4) 

and the dual integral 

/a= f /i'(z)/(z;0,J„)dz. (1-5) 

1.2 Orthogonal Rotation 

To solve the preceding expressions for h and J 2 we construct an orthonormal rotation matrix R e 
SO(n) such that the vector b is rotated to lie along the unit vector ei = (1.0, 0-). Thus, we require 
that 

Rh = |b|ei, 

where |b| is the norm of b. Thus, 



Rlt = 0i (1.6) 
V6f + • + K 

The remaining rows of R are arbitrary provided that R is an orthonormal rotation matrix. When an 
explicit construction of the remaning rows of R is needed, we construct n - 1 unit vectors ti 2 , u n 
that are mutually orthonormal and orthogonal to b. That construction can be done using Gramm- 
Shmidt or other algorithms. Nonetheless, when constructed, we set the jth row of R equal to u,. 

We now let w = J?z, and hence, z = R'w. Because R ! R = I n , it follows immediately by the linear 
transformation properties of the multivariate normal distribution that /(z; 0, I n ) is transformed into 
another JV(0, I n ) distribution /(w„;0,7„). Furthermore, the denning hyperplane of the subspace 
ft' becomes 

d - a • n - b • R'w = 0. 
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But 

h ■ R'w = (Kb) ■ w = \b\w u 
and thus, the hyperplane equation becomes 

d - a - - |b|u>i =0. (1.7) 
The transformed hyperplane now defines the new subspace 

n" = {w:ui! <(d-a-/i)/|b|}. 
The rotation also defines a new function 

h"(w) = h'(R'w) = h(TR'w + /j,). 
Letting U = (d - a • fi)/\h\, we get 

h = J _ i /(w n _ 1 ;0,J„_i)dw n _ 1 y U /i"(-w)/(wi;0,l)du> 1( (1.8) 

and similarly, 

J 2 =/" /(w n _i;0,/ n _i)dw n _i h"Mf(w 1 :0,l)dw 1 , (1.9) 
where w n _i = (w 2 , tn n ). 

2 Special Cases of ft, 

We explicitly discuss three special cases for the function h. In the first case, h is a linear function of 
x. In the second case, h is a rational function of x and in the third case h is an exponential function 
of x with either linear or quadratic terms. 

2.1 Linear Function h 

We consider the important case when the function h is a linear function of x. Thus, we can write h 
as follows 
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where q is a vector of n coefficients. Thus, we get that 



fc"(w) = h(TR'vr + /x) = c- q- /i - (JJT'q) • w. 
Let s = RT'q. Thus, Equations 1.8 and 1.9 reduce to 

h = (c - q • /i)3[-co, U] - si9 a [-oo, !7], (2.1) 

and 

J 2 = (c - q • oo] - Sl ^[U, oo], (2.2) 

where 

S .=EE^?r (2-3) 
But by Equation 1.6, R u = bj\b\, and since b = T'a and TT' = £ we get that ' 

Sl = (2.4) 

2.2 Rational Function h 

Let ft denote a polynomial function of x. As before, we can express h"(w) = h(TR'~w + ju), where 
it follows trivially that h remains a polynomial function of w under the transformation. Thus, 
it suffices to show how to evaluate h and h for the general term h" (w) = w\ 1 • ■ • w£ n where p 2 
represent arbitrary integers. We note that each ith term can be integrated separately and thus 
evaluation of 7i or h for h"(w) = w{ 1 ■ reduces to the product of terms 

where F denotes the incomplete Gamma function. 

2.3 Exponential Function h 

Let h denote the exponential function h(x) = e c_cl x . Thus, we can write 



But for any scalar s and N(Q, 1) distributed random variable x we get 
J e~" f(x;0,l)dx = e s ^ 2 J f(x;s,l)dw. 



Thus, we get that 



Thus, letting s - RT'q, and observing that s • s = q'TR'RT'q = q'Eq we get that 
I 1 = e c_q ' i+ 5 cl ' Sq Q ; [-oo,t/] 

The final case when the function h is of the form 

ft(x) = e *'* x + q 'x+c 
can also be solved, and is left as an excercise to the interested reader. 



3 The Induced Multivariate Normal Distribution 

As before, let /(x;/x, S) denote the multivariate normal density function ofx £ K° with mean 
vector fj, and covariance matrix E. Let y denote an m dimensional vector related to x by the m 
constraint equations 

J ,,-a,-x = 0, (3.1) 

where i = 1, m. In this section we derive the distribution of y induced by the distribution on x 
and the m constraint equations. 

Let p < n denote the dimension of the vector space spanned by the m vectors ai , ...,a m . We 
assume that the first p vectors a x , ...,a p are linearly independent, and therefore, that they span 
that space. Otherwise, we can reorder the vectors appropriately. In deriving the distribution of 
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y we only need to derive a distribution of y p = (yi,-~,y p ) related to x by the first p constraint 
equations in Equation 3.1. This result follows because for k > p, each y k can be expressed as a 
linear combination of the y u i < p. 

Letting A denote the p x n matrix where the ith row of A is the vector a l , we can rewrite the p 
constraint equations as 

y p - Ax = 0. (3.2) 

Let g{y p ) denote the induced distribution of y p . For a given y v , let fi(y p ) denote the subspace of 
x G 5fi n that satisfies the p constraint equations for y p . Thus, it follows that 

9(yp) = f /(x;/i,S)dx. (3.3) 

To obtain a solution for the induced distribution g we must solve the integral in Equation 3.3. 
Observe that the space fi(y p ) is the subspace generated by the intersection of the p hyperplanes de- 
fined in Equation 3.2. We precede as follows. As before, we perform a linear transformation of the 
space using the Cholesky decomposition of S, and thus transforming the distribution / (x; fx, E) 
into an jV(0, I n ) distribution. Because the transformation is linear, the space f2(y p ) is transformed 
into a new space that is the intersection of the transformed hyperplanes defined by the p trans- 
formed constraint equations. We next prove that the new space 0'(y p ) can be generated by the 
intersection of a set of p orthogonal hyperplanes. We can now rotate these orthogonal hyperplanes 
so that they lie orthogonal to coordinate axes, while leaving the Af(0,I n ) distribution invariant. 
Once rotated, we can complete the integration. 

As before, we let z = T _1 (x - p), where TT' - S denotes the Cholesky decomposition of S. 
Thus, we get p new constraint equations 

y, - a, • p - b t • z = 0, (3.4) 

where = r'a, and i = 1, ...,p. The subspace Q(y p ) is transformed into the space fi'(y p ) defined 
as the subspace of z that satisfies the constraint equations 

y p - Ap - ATz = 0. 

Furthermore, the vector z has distribution density function /(z; 0, I n ), and thus 

S(y p )=/ /(z;0,J„)dz. 

Let the vectors u! , u p form an orthonormal basis that spans the space of bi 



(3.5) 

, ...,h P . Thus, for 
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(3.6) 



where fi l3 = h t ■ Mj . Let A denote the p x p matrix of j3,j. Let U denote the p x n matrix with ith row 
containing u { for all i = 1, Thus, we can rewrite Equation 3.4 as 

y p -Afi- AUz = 0, (3.7) 

or, since A is invertible, as 

A _1 y P - A _1 A^ - Uz = 0. (3.8) 

Now let w - Rz where R is an orthonormal rotation such that Rut = e, for i = 1, ...,p where e z 
denotes the ith unit basis vector. Observe that UR' is the p x n matrix with ijth element equal to 
S v for j <i, and zero otherwise. Thus, UR'-w = w p where w p = ...,w p ). Substituting z = i?.'w 
in Equation 3.8 we get 

w p = A- 1 y P -A- 1 A/i. (3.9) 

Thus, the subspace fi' (y p ) is transformed into the subspace 

n"(y P ) = {w : w p = A-^p - A' 1 Ap,}, 

and the density function /(z; 0, J„) is transformed into the density function /(w; 0, /„). Thus, we 
get that 

g(y P ) = J /(w ; o,J n )dw 
= /(y p ;AA',^) : 

where the last equality follows from well known properties of the multivariate normal density func- 
tion. Thus, the density function for y p induced by the multivariate normal distribution on x and 
the p independent constraint equations is another multivariate normal distribution with covariance 
matrix AA' and mean vector A\i. 



4 Expectation of General Functions Over Polyhedral Subspaces 



We revisit the integral h in Equation 1.1. We are now interested in the solution of h when the 
subspace of integration ft is a polyhedra described as the region enclosed by the intersection of 
multiple halfspaces, 

di-a.-xX), (4.1) 

for i = 1, ...,m, or equivalently, 

d - Ax. > 0, (4-2) 

where A is an m x n matrix with ith row equal to a,. Let O t = {x : d l - a t • x > 0}, and therefore, 
O = ! 0,. Observe that we allow the d, and a 2 to be negative, that is, 

-|d I | + |a I j-x>0, 

and hence, we may include the complement halfspaces ft, = {x : \d,\ - |a,| ■ x < 0} in the definition 
of the polyhedra ft. (Observe that in the preceding expressions, |a,,| represents the vector a d with 
all positive entries, and not its norm.) 

Proceeding similarly to the derivation in Section 1, we perform a Cholesky decomposition TT = 
S, enabling us to reduce Ji to 

h = ( h'(z)f(7.;0,I n )dz (4.3) 
J w 

where h'{%) = h{Tz + n) and 0' = n™ where 

Q[ = {% : d i - a, • p - h % ■ z > 0}, 

andb, = T'a,. 

Equation 4.3 cannot be solved exactly in closed form. In fact, we show in Section 6 that the solu- 
tion of h is NP-hard. We can approximate h as follows. We identify a set of mutually orthogonal 
halfspaces that best approximates the polyhedra fl'. Let 

fi'/ = {z : a - u x • z > 0}, 

for i = 1, ...,r, denote those r orthogonal halfspaces, and let fi" = n[ =i nj' denote the polyhedra 
defined by the orthogonal halfspaces. By construction, the vectors u t are orthogonal. Furthermore, 
we assume that they have been normalized, and hence are orthonormal vectors. We can express ft" 
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as the solution of 



c - Uz > 0., (4.4) 

where V is the r x n matrix with ith row equal to vl z . Observe that r < p, where p is the dimension 
of the subspace spanned by the vectors h z (or equivalently, by the vectors a,}. With the preceding 
approximation of the polyhedra SI' , we have that 

I[= f /i'(a)/(z;0,/„)<fa (4.5) 

Jq" 

approximates h • We now show how to solve I{ exactly. In the next section we show how to 
construct the halfspaces f2". 

We construct a rotation matrix R e SO(n) such that i?u, = e 4 for i = 1, r where e 4 is the ith 
unit vector. For i = 1, r, the ith row of R is equal to u ; . For i > r, the rows of R are arbitrary 
subject to the constraint that J? is a rotation matrix. When explicit construction of those rows is 
necessary, we construct n - r orthonormal vectors v, that are orthogonal to the vectors u,. That 
construction can be done using the Gram-Schmidt algorithm. The ith row of R for i > r is set equal 
to the vector Vj_ r completing the construction of R. 

Let w = R'z. This rotation transforms the constraints in Equation 4.4 into 

c - UR'w > 0. (4.6) 

But by construction of R, UR' is the r x n matrix with ijth element equal to <5 £j for j < i and zero 
otherwise. Thus, the halfspaces fJ' 2 ' are transformed into fl 2 = {w : w, < c 2 }. Furthermore, the 
density function /(z;0,/ n ) is transformed into /(w;0,J n ) and the function h' becomes h"{w) = 
h'(R'w) = h{TR'w + n). Thus, I[ becomes 

^ = J mn r f[ /" ^"( w )/( w ; °. J ») dw (4J) 

4 . 1 Construction of Halfspaces Q." 

Recall that we desire to construct a set of orthogonal halfspaces O^', % = 1, r, that best approxi- 
mate the polyhedra W enclosed by the intersection of the m halfspaces 

U[ = {z : d, - a, - jj, - b, • z > 0} 

Let v 1 , i - 1, ...,p denote an orthonormal vector basis that spans the p-dimensional subspace de- 
fined by the vectors h t ,i = 1, m. Let b\ denote the normalized vectors h t . We can express each 
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b, as 

where A y = b t • v., . To obtain a good approximation of the polyhedra denned by the ft; using 
the orthogonal halfspaces ft'/ we would like the coefficients \ t} to be either large (close to one) 
or small (close to zero). For example, if for each i, X tj equals one for some j and zero otherwise, 
then the orthogonal halfspaces ft/ exactly reproduce the polyhedra ft'. (This result trivally follows 
by noting that the b, already form an orthonormal set of vectors and hence the halfspaces ft' ; are 
orthogonal.) To achieve an optimal set of coefficients A„ we can rotate the vectors v 3 so that they 
are approximately oriented along the bj. 

Let A denote the m x p matrix with ijth element equal to \, . We rotate the vectors vj by a matrix 
R such that the new coefficients A* = AR are either small or large. The rotation can be found by 
minimizing the well known criterion 

V = tpK) 4 -it\P^. ' '<4-8) 

Equation 4.8 is well known in factor analysis, and minimization of V when 7 = 1 yields the varimax 
rotation and when 7 = 0 yields qartimax rotation. 

We now define the new orthonormal basis = R'v t , i = l,...,p, and thus 

h^jz^- (4 - 9) 

We approximate each b, by the Uj with largest coefficient A^ . In general, if several b\ are approxi- 
mately collinear, then we do not need the entire basis u ; , j = 1, p to approximate the vectors b t . 
Assume that r <p basis vectors are required. We further assume that the vectors u 3 have been re- 
ordered so that Ul , . . . , u r define the desired r basis vectors. We define the mapping a from the index 
set {1, ...,m) into the index set {1, ...,r} such that u ctW is the basis vector that best approximates 
bj. We can now define the m orthogonal halfspaces ft" as 

ft;' = {z -. d, - a, • fi - |b ; |u CT(;) • z > 0} (4.10) 

Note that the goodness of this approximation rest largely on how well the varimax or quartimax 
rotations can generate a matrix A* with entries that are either close to one or close to zero. 

Because r < m, all vectors b, such that i e a' 1 (j) will get mapped into the same basis vec- 
tor u,. Thus, those halfspaces ft" will have parallel boundaries, and it follows trivially that their 
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intersection is the halfspace where k e a 1 (j) is the index given by 
argmin l6(T - 1(j) . 

Thus, for Equation 4.4, we get 

C, = mm rr-j 

and U is the r x n matrix with ith row equal to uj. 

5 Special Cases of h for Polyhedral Integration 

We explicitly discuss the cases when h is a linear function of x and when it is an exponential func- 
tion of x. 

5.1 Linear Function h 

As before, we express h as 

Mx)=/-q-x 

where q is a vector of n coefficients. 

7»"(w) = h(TR' W + p) = /- q- A* - CRT'q) - w. 
Let s = ET'q. Thus, Equation 4,7 reduces to 

^=(/-q-M)nO[-~, Cl ]-E Sj ^[-oo )Cj ] f[ ( 5 -D 

2=1 J=l 2=1. 

where for j < r 

^EE^^^ (5 - 2) 

since by construction , the j th row of .R is u, for j < r. Thus, if we let s p denote the vector containing 
the first p components of s, we get 

s p = E/T'q. 
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(4.11) 



(4.12) 



In the special case when the vectors h Jt j = I, ...,m (and hence r = m) are mutually orthogonal, 

we get 




(5.3) 



since hj = T'a. 3 and TT' — E. 

5.2 Exponential Function h 

Let h denote the exponential function fe(x) = e- f ~ q x . Thus, we can write 

A"(w) = J-** ■ e-( RT '«>™. 
Thus, letting s = i?T'q, and using Equations 2.6 and 2.7 we get 

I[ = e'~* ^ + i q ' Sq f[ 9[-oo, Ci] - • (5.4) 

where we used the relationship s • s = q'Ti?'j?X"q = q'Sq. 

6 Integration of the Multivariate Normal Distribution over Poly- 
hedral Subspaces is NP-hard 

Much research has been dedicated to solve, or efficiently approximate, the integral /(x; /i, E) over 
an arbitray region. Even the simple bounding region x t < 0 for i = 1, ...,n, has eluded mathe- 
maticians. We prove that the problem of integrating even /(x; 0, I n ) over an arbitrary polyhedra is 
NP-hard. The proof consists of a straightforward reduction of the NP-hard problem of computing 
the volume of a polyhedra. Let a, - x < 0, % = 1, ...,p define the polyhedra O. The volume enclosed 
by that polyhedra, vol(ft) can be easily shown to equal 

vol(fi) = / /(x;0,/ n )dx, 
Jn 

thus proving our claim. 



12 



7 Appendix 1 



We present some well known results on the integration of a normal distribution. Let f{x\/j.,a 2 ) 
denote the normal distribution with mean p. and variance a 2 . Let q = q{x) = (x - ii)ja. Then, 



j f(x:p,a 2 )d X ^^Eri |j=] , 



We let 9f[u.v] and denote, respectively, the integrals in Equations 7.2 and 7.1 evaluated 

between the limits u < x < v. Also observe that these integrals trivially satisfy the following 
properties 

9[ U> «] = 3[-oo,i;]-9[-oo ) u] 
3[-oo,*] + 3[ U ,oo]=l 
9 as [«,«] = a a! [-oo,«]-9 s [-oo,«] 
9 x [-oo,u] + 3[u, oo]=/i 

We now evaluate the derivatives of the integrals 9 and $s x . From the above properties, it suffices to 
evaluate the integrals Q[— oo, v] and oo, v]. Letting q = q(v) = i ^ L , 
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1 Introduction 



This report focuses on several incarnations of Rapt Technologies proprietary uniform pro- 
duction policy for production planning and component allocation under scarce resources. In 
this report we revisit the definition of the policy for the MCMP model (previously defined 
in [1,2]) and we derive closed form analytical solutions to the expectation of a multivariate 
normal distribution under this policy. What we are seeking in this paper is an expression 
for the expected production vector given a demand distribution (which in general, will not 
be fully satisfied in all situations). 

2 Background 

Let /'(x; jj., S) denote the multivariate normal density function of x e 3ft" with mean vector 
(i and covariance matrix E. Let [m] = {1, rn}. Let h{x) denote an arbitrary multivariate 
functions of x and let Q. denote the polyhedra described as the region enclosed by the 
intersection of multiple half spaces, 

d l -a i -x>0, (2.1) 

for i G [to]. Equivalently, define 

d - Ax > 0, (2.2) 

where A is an m x n matrix with ith row equal to a 2 . Let Q t = {x : d % — a,. - x > 0}, and 
therefore, fi = n t£ [ m ]fi,. 

The space Q, denotes the feasible region where the component level d is sufficient to 
meet any demand x in that space. When x falls outside of that space, then we must in- 
voke either a production policy that prescribes the level of products to manufacture while 
satisfying the component levels, or equivalently, an allocation policy that prescribes how to 
redistribute the components among the products. The uniform policy represents a produc- 
tion policy, rather than an allocation policy. In general, production policies are preferred 
over allocation policies since company revenue is linked directly to product manufacture 
and not component consumption. 
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A production policy y is a mapping from the space f2 into the space Q. such that for any 
x € Ct we have that y(x) < x. A maximal production policy is a production policy that maps 
O into the boundary subspace of Q. Thus, a maximal production policy will deplete at 
least one component. 

InTR-02-99 we obtain closed-from expressions for the following general integral 

1(d) = f x/(x;/*,E)dx, (2.3) 
J SI 

and 

7(d) = / x/(x;/x,£)dx (2.4) 
J si 

by noting that 

1(d) = Ex -1(d) = a* -1(d). " (2.5) 

TR-02-99 also solved the above integrals for the priority allocation scheme (or equivalently, 
the priority production policy). In this report we extend the development to include other 
production policies. 

2.1 Using these results 

While this paper only presents the calculation of the expected production vector, this result 
can be directly used to evaluate the expected value of any linear function of the production. 

3 The u-policy 

Let u denote some arbitrary normalized vector. The (u-policy) y is defined for all as 

ice 



Jx-a(x)u ifxe 

y(x) = < 

[ x otherwise 
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for some scalar valued function a such that 



d - Aq, > 0, 



(3.2) 



and 



q(x) < x. 



(3.3) 



The choice of vector u and function a are arbitrary, but in general, they can be chosen to 
optimize some objective function. 

The most important characteristic of the w-policy is that it prescribes a particular di- 
rection for the x - q(x) vector, irrespective of the value of x. This can lead to unexpected 
prescriptions, such as a negative production corresponding to a positive demand vector. 

3.1 Transformation of the Integral i 

In TR-02-99 we show that after several transformations, the general integral I(h o q, d) can 
be reduced to the following equivalent integral 




(3.4) 



where h"(w) = h(q(T J R'w + Af))> 




(3.5) 



and a is the mapping from the index set [m] into the set [r] denned in TR-02-99. 



3.2 Determining of a 



Before we show how to compute the integral 




(3.6) 
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we need to define a. For linear, polynomial, or exponential forms of a this computation 
can be done using the methods derived in TR-02-99. For the uniform production policy 
however we choose 

a(x) = max [ &i X — 1 - (3-7) 

,e[m] [ a, • u J 

This choice guarantees that for any x e Cl and a given u, the production level q(x) is the 
highest level compatible with the restrictions in Equations 3.2 and 3.3. 

Using the transformations of TR-02-99 we can reduce or(x) as follows. First, letting 
x = Tz + n we get 

where b, = T'a 2 . Next, from the orthomax rotation used to 'orthonormalize' the spape Q. 
we have that the m vectors h z can be approximated by the r orthonormal vectors u 3 such 
that bj ~ |bi|u ff(i ). As before, a is a mapping from the index set [m] into [r]. Substituting 
for the b, and defining the rotation matrix R such that Ru a{l) = e a(i) , where are unit 
vectors, we get 

a"(w) = a(#w) = max [Au VW + li] (3-9) 
where, noting that jb t | — y / a^£a;, 



We can rewrite a" a 



3.3 Computing / 

We can now show how to compute 

1= f a(x)/(x;/*,E)dx. (3.11) 

Jn. 

We use the results of TR-02-99 and of the preceding section to reduce I to 

T = f[ f C ' a"(w)/(w;0,J u ;kiw, (3.12) 

where 

_ . \ dj - a, ■ /i l 

Cz ,™ ? WL |b,| J 

as denned in TR-02-99. It now follows directly from the results in the Appendix that 

I = Y^r max (0 jWi + lo )Y[^[-oo,xnm(c k ,Uik(x))}f(x-,OA)dx (3.13) 
where, letting 

Mj = max WjWt + jj), 

we have that 

\Mi -7,1 
17,* = min ' fJ . 
jeo-M*) L & J 

The integral I is a sum of r < m univariate integrals. The integral can be readily solved 
using numerical integration techniques such as quadratures or even Simpson approxima- 
tions. Note that the integrand is a product of functions 3f[] that results in a well behaved 
monotonic increasing function. 
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4 The x-policy 



The x-policy is similar conceptually to the u-policy, except that demand is scaled propor- 
tionally: 

^ fa(x)x ifxeO (41) 

!x otherwise 

The x-policy is superior to the u-policy that it ensures that positive demand maps to 
positive production. 

4.1 Determining of a for the x-policy 

Alpha is chosen similarly to the u-policy case: 

a(x) = min f-A_l . (4.2) 

As before, first, we let x = Tz + jj,: 

■^-"^"'-ffitbTiTI^] (4 ' 3) 

where b l = T'a,- Next, from the orthomax rotation used to 'orthonormalize' the space Q 
we have that the m vectors b 4 can be approximated by the r orthonormal vectors xx 3 such 
that hi ~ |b J |u (T( , ) . As before, a is a mapping from the index set [m] into [r]. Substituting 
for the bj and defining the rotation matrix R such that Ru a ^ = e^), where e, are unit 
vectors, we get 

a"(w) = «'(ii'w) = min | — 1 (4-4) 

y J te[m] lfaw a{{) +7tJ 

where, noting that |b t | = ^/a^Eai, 

= V^pl (4.5) 
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and 



(4.6) 



We can rewrite a" as 



a"(w) - min [ min (— — t ) 1 . (4.7) 

k ; ie[r] bea-Hxy P 3 w t + 7/ J 



4.2 Evaluating the expected production vector for the x-policy 

In the x-policy case, the integrand that is of concern is a(x)x. A for the u-policy, we attempt 
to reduce the multidimensional integral to a sum of one-dimensional ones. However, we 
need to make the following approximation: 

o/'(w) = min [- - ] * ^ — V '(4-8) 

The approximation here is exact when Piiv a(l ) + 7i > 0 for a11 * e H- However, it is 
potentially very far off otherwise. We next show that the approximation holds in the region 
of importance for the integral I. 

Assume that the region of appreciable demand probability falls completely in the 
quadrant Q = {x|a;j > CM}, and that all products require a nonnegative amount of each 
component {a^ > (M, j). In this case, we have that 

ftti) ffW + 7l = a, ■ x > OVz, (4.9) 

and hence the approximation in Equation 4.8 holds for the purposes of evaluating the 
integral I. 

With the above transformation, the result in the Appendix applies, and the integral 
can be evaluated as follows: 



= f a(x)x/(x;/i,S)dx. (4.10) 
Jn 
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We use the results of TR-02-99 and of the preceding section to reduce I to 

I = Cr dw/(w;0,/ n )a"(w)(T J R'w + Ax) (4.11) 

= Ti?'Ii+J 2 /i, (4-12) 
In. = TR' p Cr dw/(w;0 1 7„)Q"(w)w, (4.13) 

h = Cr ^w/(w;0,I„)a"(w) (4-14) 

(4.15) 



« = min • (4-16) 

The value of I 2 now follows directly from the results in the Appendix: 

J 2 = tr fMpdw Y[Q[-<x>,"**{ckMw))] (4.17) 

where, letting 

M, = max (/9j'to 2 +7j), 

we have that 

The value of Ii, which is an r-dimensional vector, can be obtained with a slight modi- 
fication of the result in the Appendix. Using the same notation as above, for all I € [r] 

if) = y IT fMVdw 9g[ ^ n(ct)l7 ., (ui))] Yl csl-oo.mMckMv 

^ 1 ./-co + lj) k % 



max jetr -i (0 (£,u> +7j) *i 



The integrals 1^ and J 2 are sums of r < m univariate integrals, which can be read- 
ily solved using numerical integration techniques such as quadratures or even Simpson 
approximations . 

5 The Local w-Policy 

The u- and x-policies prescribe a demand-production mapping that is manifestly irrational 
in some cases. In particular, the problem of "intercomponent" effects arises where a com- 
ponent that gates production for one product will diminish production over all. This arises 
because the previous policies compute a single a for all products. 

On an event by event bases, an "iterated policy" represent rational demand-production 
mappings, but analytic formulas are not available for the expected production over these 
policies. A compromise is to compute a separate a for each product. This we refer to as 
the local u -Policy. 

As before, for product i, let V(i) denote the components in the BOM of i. Also, as 
before, for any component j we can define the subspaces Qj = {x : dj - a 3 • x > 0}. We 
now define the subspaces O l for each product i as 

If x e O. 1 then component allocation is sufficient to meet product demand for the ith prod- 
uct, even though it may not be sufficient to meet demand for other products. 

We make the important observation that in computing expected values for functions 
h(x) that can be separated into a sum, 

fc(x) = h l {x l ) + --- + h n {x r) ) 

the appropriate subspace of integration for the feasible production for each function h 3 is 
9J. 

For any subspace Q, let En denote the expectation over support fi. Thus, the expecta- 
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tion of h can be expressed as 

Eh = 1t [E^Ms.) + E nt ft,{ 9l (x)] 

where q(x) is any production policy that maps each Xi £ Cl l into <j;(x) e 
We define the local zt-policy q as 



.w = { ; 



if x e sv, 

- a t (x, d)u 4 otherwise, 



a z (x, d) = max 



j6X>(0 [ u J 

We now show that with the preceding definition, for any k, we have that d k — & k -q > 0. 
Substituting q we get 

d k - & k - q = d k - a fc ■ x + ]P a kl u t i 



jer>(i) I a.j • u J 
But observe that a ki — 0 if k g T>(i), and therefore 

\aj ■ x - ciJ a*; • x- d fe 

max — > , 

3 €V(i) [ SLj U J a k ■ u 

proving that d k — a^ • q > 0. 

6 The orthomax policy 

The orthomax production policy prescribes that for any infeasible demand x e Q the pro- 
duction should be given by transforming x to w-space, finding the nearest point on the 
production boundary, and transforming the boundary point back to x space. Using the 
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notation of the previous sections, this algorithm translates to: 

fri?'w'(w(x)) + (j, + U(x) where x e fi 
q(x) = ^ (6.1) 
! x otherwise, 

where w(x) = i?T _1 (x - //) and 

w'(w) = ^min(c 2; iu 2 )e I = (min(ci, u>i), min(c 2 , w 2 ), - ■ ■ , min(c r , w r )) . (6.2) 

U(x) is the affine map given by the equation 

U(x) = n + TR'RT- 1 ^ - ft). (6.3) 

Conceptually, U(x) is the affine map that gives the component of x mapped to 0 by the 
affine map x i-> w. 

The expected production under the orthomax policy is given by the integral 



f dx/(x;^,S)q(x) 

f dxx/(x;^,S)+ /<fx/(x;/t,E)(q(x)-x) 
inun Jn 



where Ii is the integral over Q. on the second line of the formula above. 

In order to make progress evaluating 1\, we change the region of integration Q. to the 
complement of the "orthomax region" given by Q ! = {x| (i?T _1 (x — p))-e,; < c* V? G [r]}. 
At the same time, we redefine q(x) to be the equal to the non-identity map on ft', which is 
another small, but appreciable approximation. With these approximations, we can write 

I 1= /dx/(x;^,5])(q(x)-x)« / dx/{x; /*, E)(q(x) - x). (6.4) 



n 



A closed form solution can be now obtained: 

Ii = TR ' J""' Cr dw/(w;0,/„)(w-c) + j _ dx/(x; /x, E)U(x) (6.5) 
= ^'((^(c^^-d^c^cx))),...,^^,^)-^^,^))) (6.6) 

(6.7) 

where we made use of the fact that since TJ(^t + £) = — U(/x — £) (i.e U is odd around 
x = /n), /(x; /it, S) is even around x = /u., and U(x) • a, = OVi, we have 

jf dx/(x;/i,E)U(x) = 0. (6.8) 



7 Appendix 

We show how to reduce the following multivariate integral into a sum of univariate in- 
tegrals that can be readily computed. Let ce(x) = G(maoc ie [ n j + 7,]), where G(x) is a 
well-behaved 3£ — >• 3F£ function. We desire to compute the integral 

1 = (fl /_^) (x; 0, J n )a(x) . (7.1) 

Consider the case when fi\x\ + 71 is maximum. This occurs when for all j 7^ 1 /^cc,- + 7^ < 
/3ix x + 71, or when 

A^l + 71 ~ Ij 
Xj ~ Pi 

Define 
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Thus we get that 

I = ^f\ 1 /( Ii :0,l)G(A Il l7 l )nr (C, '' ,,,W) /( 1 j; 0 . 1 )^ 
= jr r tef{x#A)G{&x + li )Jl%[-<x>,T^{c ] ,U t3 {x))] 
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